R version 4.0.3 (2020-10-10) – “Bunny-Wunnies Freak Out”
Packages used for NMDS: vegan (version 2.5-7)
The document serves as an example of analyses that will be conducted to identify natural variations in benthic communities across Virginia. These NMDS will support the Genus level IBI development process. This analysis is the first run of all of reference sites in Virginia. No West Virginia DEP data is used in this analysis. Several reference sites have been sent back to Virginia DEQ biologists for one more review.
The dataset used is a subset of reference stations collected in Virginia. If stations appeared in the dataset more than 4 times, then the most recent 4 samples were used and the rest removed. Taxa that occurred in the dataset < 3% of the time were removed. The data was log10 +1 transformed. Environmental factors were compiled for each station and used to plot over the NMDS to show environmental variation associated with the community matrix. The envfit function in Vegan was used to plot the continuous environmental variables.
The first step was to read in the reference site bug taxa list and environmental factors dataset for each station. Join the environmental dataset with the bug dataset to account for multiple observations of each station and collection date and time.
Check to make sure the bug and environmental join was successful:
Number of rows in Community Matrix: 854
Number or rows in Environmental Matrix: 858
The data was log10+1 transformed. Rare taxa (<=3%) were removed.
## Run 0 stress 0.1677427
## Run 1 stress 0.1679284
## ... Procrustes: rmse 0.004018228 max resid 0.07279403
## Run 2 stress 0.1678503
## ... Procrustes: rmse 0.002571914 max resid 0.05503484
## Run 3 stress 0.1679616
## ... Procrustes: rmse 0.003128083 max resid 0.06762832
## Run 4 stress 0.1680443
## ... Procrustes: rmse 0.003386152 max resid 0.06767626
## Run 5 stress 0.1677157
## ... New best solution
## ... Procrustes: rmse 0.002936922 max resid 0.07307232
## Run 6 stress 0.168044
## ... Procrustes: rmse 0.004257538 max resid 0.07322163
## Run 7 stress 0.1679439
## ... Procrustes: rmse 0.002916938 max resid 0.06766017
## Run 8 stress 0.1677183
## ... Procrustes: rmse 0.0004981453 max resid 0.0125922
## Run 9 stress 0.1677723
## ... Procrustes: rmse 0.001892368 max resid 0.045066
## Run 10 stress 0.1680416
## ... Procrustes: rmse 0.003193013 max resid 0.06764242
## Run 11 stress 0.1677989
## ... Procrustes: rmse 0.001448638 max resid 0.04201973
## Run 12 stress 0.1679346
## ... Procrustes: rmse 0.003438932 max resid 0.06685484
## Run 13 stress 0.1678304
## ... Procrustes: rmse 0.002091694 max resid 0.05492619
## Run 14 stress 0.1677237
## ... Procrustes: rmse 0.002083582 max resid 0.04472721
## Run 15 stress 0.1689301
## Run 16 stress 0.1688288
## Run 17 stress 0.1677408
## ... Procrustes: rmse 0.002531716 max resid 0.07312253
## Run 18 stress 0.1678179
## ... Procrustes: rmse 0.002617366 max resid 0.04820438
## Run 19 stress 0.1677173
## ... Procrustes: rmse 0.001751015 max resid 0.04508061
## Run 20 stress 0.168979
## Run 21 stress 0.1679634
## ... Procrustes: rmse 0.002784021 max resid 0.06755936
## Run 22 stress 0.1679618
## ... Procrustes: rmse 0.003705789 max resid 0.0730803
## Run 23 stress 0.1677402
## ... Procrustes: rmse 0.002054085 max resid 0.04526348
## Run 24 stress 0.167717
## ... Procrustes: rmse 0.0003823337 max resid 0.009819117
## ... Similar to previous best
## *** Solution reached
##
## Call:
## metaMDS(comm = NMDSone[, 6:213], k = 3, trymax = 1000)
##
## global Multidimensional Scaling using monoMDS
##
## Data: NMDSone[, 6:213]
## Distance: bray
##
## Dimensions: 3
## Stress: 0.1677157
## Stress type 1, weak ties
## Two convergent solutions found after 24 tries
## Scaling: centring, PC rotation, halfchange scaling
## Species: expanded scores based on 'NMDSone[, 6:213]'
Envfit results from Vegan package :
##
## ***VECTORS
##
## NMDS1 NMDS2 r2 Pr(>r)
## JulianDate -0.76892 -0.63934 0.0738 0.656
## Latitude 0.93367 0.35813 0.2793 0.170
## Longitude 0.88025 -0.47450 0.2551 0.201
## US_L3CODE 0.40994 0.91211 0.3416 0.105
## Order -0.84033 0.54207 0.7786 0.002 **
## EDASOrder -0.84033 0.54207 0.7786 0.002 **
## totalArea_sqMile -0.42748 0.90403 0.4581 0.048 *
## ELEVMEAN -0.63475 0.77272 0.5463 0.012 *
## SLPMEAN -0.11679 0.99316 0.7957 0.002 **
## wshdRain_mmyr -0.40410 0.91472 0.4617 0.044 *
## siteRain_mmyr -0.73238 -0.68090 0.2914 0.134
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 999
##
## ***FACTORS:
##
## Centroids:
## NMDS1 NMDS2
## SeasonFall -0.7303 -0.5855
## SeasonOutside Sample Window -0.5972 -0.0876
## SeasonSpring -0.6775 -0.2603
## GradientMACS -0.6944 -0.3874
## CoastalTRUE -0.6944 -0.3874
## US_L3NAMEMiddle Atlantic Coastal Plain -0.6942 -0.7009
## US_L3NAMESoutheastern Plains -0.6945 -0.2132
## US_L4NAMEChesapeake-Pamlico Lowlands and Tidal Marshes -0.6963 -0.3463
## US_L4NAMEMid-Atlantic Flatwoods -0.6929 -0.9374
## US_L4NAMERolling Coastal Plain -0.6945 -0.2132
## ASSESS_REGPRO -0.6949 -0.2374
## ASSESS_REGTRO -0.6929 -0.9374
## TidalN -0.6941 -0.3942
## TidalT -0.6963 -0.3463
## VAHUSBCB -0.6220 -0.5269
## VAHUSBCU -0.7650 -0.5287
## VAHUSBYO -0.6609 -0.0358
## BasinChes. Bay and Small Coastal Basin -0.6220 -0.5269
## BasinChowan and Dismal Swamp River Basin -0.7650 -0.5287
## BasinYork River Basin -0.6609 -0.0358
## Basin_CodeChowan-Dismal -0.7650 -0.5287
## Basin_CodeSmall Coastal -0.6220 -0.5269
## Basin_CodeYork -0.6609 -0.0358
## CountyCityNameCaroline -0.7402 0.0413
## CountyCityNameDinwiddie -0.8776 -0.3233
## CountyCityNameHanover -0.5817 -0.1129
## CountyCityNameKing and Queen -0.5477 -0.7076
## CountyCityNameNorthumberland -0.6963 -0.3463
## CountyCityNameSouthampton -0.6929 -0.9374
## CountyCityNameSussex -0.7565 0.2866
## WQS_CLASSIII -0.7259 -0.0647
## WQS_CLASSVII -0.6769 -0.5666
## WQS_SPSTDS -0.6753 -0.3451
## WQS_SPSTDSNEW-21 -0.9431 -0.9367
## WQS_PWS -0.6944 -0.3874
## WQS_TROUT -0.6944 -0.3874
## WQS_TIER_III -0.6944 -0.3874
##
## Goodness of fit:
## r2 Pr(>r)
## Season 0.1787 0.386
## Gradient 0.0000 1.000
## Coastal 0.0000 1.000
## US_L3NAME 0.2986 0.043 *
## US_L4NAME 0.4624 0.027 *
## ASSESS_REG 0.4510 0.004 **
## Tidal 0.0015 1.000
## VAHUSB 0.2919 0.131
## Basin 0.2919 0.131
## Basin_Code 0.2919 0.131
## CountyCityName 0.8748 0.004 **
## WQS_CLASS 0.3192 0.036 *
## WQS_SPSTDS 0.1529 0.282
## WQS_PWS 0.0000 1.000
## WQS_TROUT 0.0000 1.000
## WQS_TIER_III 0.0000 1.000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 999
##
## 840 observations deleted due to missingness
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$Season, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## Fall Outside Sample Window Spring
## delta 0.3535 0.3975 0.4537
## n 414 13 427
##
## Chance corrected within-group agreement A: 0.02743
## Based on observed delta 0.4043 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$US_L3NAME, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## Blue Ridge Central Appalachians Middle Atlantic Coastal Plain
## delta 0.2731 0.3524 0.5705
## n 154 40 25
## Northern Piedmont Piedmont Ridge and Valley Southeastern Plains
## delta 0.352 0.4253 0.3165 0.5478
## n 101 138 259 137
##
## Chance corrected within-group agreement A: 0.0938
## Based on observed delta 0.3767 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$Basin_Code, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## Appomattox Chowan-Dismal James-Lower James-Middle James-Upper New
## delta 0.4055 0.5345 0.5585 0.3775 0.3011 0.3248
## n 15 33 36 61 120 78
## Potomac-Lower Potomac-Shenandoah Rappahannock Roanoke Small Coastal
## delta 0.4911 0.337 0.3983 0.4073 0.513
## n 31 42 144 87 21
## Tennessee-Big Sandy Tennessee-Clinch Tennessee-Holston Yadkin York
## delta 0.3435 0.3075 0.2467 0.1743 0.5018
## n 16 50 53 4 63
##
## Chance corrected within-group agreement A: 0.07742
## Based on observed delta 0.3835 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$ASSESS_REG, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## BRRO NRO PRO SWRO TRO VRO
## delta 0.3482 0.3944 0.54 0.3084 0.5523 0.3624
## n 250 195 122 156 21 110
##
## Chance corrected within-group agreement A: 0.07205
## Based on observed delta 0.3857 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$Gradient, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## MACS Riffle
## delta 0.5516 0.3429
## n 163 691
##
## Chance corrected within-group agreement A: 0.07933
## Based on observed delta 0.3827 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$Order, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## 1 2 3 4 5 6
## delta 0.4131 0.4177 0.4227 0.423 0.3213 0.1486
## n 277 248 172 114 41 2
##
## Chance corrected within-group agreement A: 0.007247
## Based on observed delta 0.4127 and expected delta 0.4157
##
## Significance of delta: 0.003
## Permutation: free
## Number of permutations: 999
##
## Call:
## mrpp(dat = bugsnms[, 6:214], grouping = samplescoresenv$WQS_CLASS, distance = "bray")
##
## Dissimilarity index: bray
## Weights for groups: n
##
## Class means and counts:
##
## III IV V VI VII
## delta 0.4941 0.3405 0.3399 0.2468 0.5416
## n 331 185 87 209 42
##
## Chance corrected within-group agreement A: 0.06913
## Based on observed delta 0.387 and expected delta 0.4157
##
## Significance of delta: 0.001
## Permutation: free
## Number of permutations: 999